/**
 * 
 */
package ua.lviv.franko.solvers.stressanalysis.oneDimention.cylindrical;

import java.util.ArrayList;
import java.util.Iterator;
import java.util.List;

import org.jfree.data.xy.XYSeries;
import org.jfree.data.xy.XYSeriesCollection;

import ua.lviv.franko.integration.ConstantFunction1D;
import ua.lviv.franko.integration.FunctionWithoutDerivatives1D;
import ua.lviv.franko.integration.GaussQuadrature;
import ua.lviv.franko.integration.IFunction1D;
import ua.lviv.franko.matrix.RowType;
import ua.lviv.franko.matrix.hash.lu.HashMatrixLUStyle;
import ua.lviv.franko.matrix.hash.lu.LUFactorization;
import ua.lviv.franko.visualisation.NewOkCancelDialog;

/**
 * @author marko
 * 
 */
public class CylindricMembraneSolver {
	public static final int					NUMBER_OF_UNKNOWNS	= AbstractElement1D.NUMBER_OF_UNKNOWNS;
	public static final double				EPS					= 0.000000000001;

	protected ArrayList<AbstractElement1D>	elementsList;
	protected HashMatrixLUStyle					mainMatrix;
	protected double[]						b;
	protected double[]						result;

	public CylindricMembraneSolver(ArrayList<AbstractElement1D> elementsList) {
		this.elementsList = elementsList;
	}

	public void generate() {
		List<Integer> lastElemetNumbers = elementsList.get(elementsList.size() - 1).getNodesNumbers();
		int N = lastElemetNumbers.get(lastElemetNumbers.size() - 1) + 1;
		mainMatrix = new HashMatrixLUStyle(RowType.ARRAY_TYPE, N * NUMBER_OF_UNKNOWNS);
		b = new double[N * NUMBER_OF_UNKNOWNS];

		for (Iterator<AbstractElement1D> iterator = elementsList.iterator(); iterator.hasNext();) {
			AbstractElement1D element = iterator.next();
			// element.generateAll();
			ElementMatrix data = element.generateAll();
			List<Integer> pos = data.getNodeList();
			for (int i = 0; i < pos.size(); i++) {
				for (int j = 0; j < pos.size(); j++) {
					int pos1 = NUMBER_OF_UNKNOWNS * pos.get(i);
					int pos2 = NUMBER_OF_UNKNOWNS * pos.get(j);
					mainMatrix.setValue(pos1, pos2, mainMatrix.getValue(pos1, pos2) + data.matrix[NUMBER_OF_UNKNOWNS * i][NUMBER_OF_UNKNOWNS * j]);
					pos2++;
					mainMatrix.setValue(pos1, pos2, mainMatrix.getValue(pos1, pos2) + data.matrix[NUMBER_OF_UNKNOWNS * i][NUMBER_OF_UNKNOWNS * j + 1]);

					pos1 = NUMBER_OF_UNKNOWNS * pos.get(i) + 1;
					pos2 = NUMBER_OF_UNKNOWNS * pos.get(j);
					mainMatrix.setValue(pos1, pos2, mainMatrix.getValue(pos1, pos2) + data.matrix[NUMBER_OF_UNKNOWNS * i + 1][NUMBER_OF_UNKNOWNS * j]);
					pos2++;
					mainMatrix.setValue(pos1, pos2, mainMatrix.getValue(pos1, pos2) + data.matrix[NUMBER_OF_UNKNOWNS * i + 1][NUMBER_OF_UNKNOWNS * j + 1]);

				}
				b[pos.get(i) * NUMBER_OF_UNKNOWNS] += data.vector[NUMBER_OF_UNKNOWNS * i];
				b[pos.get(i) * NUMBER_OF_UNKNOWNS + 1] += data.vector[NUMBER_OF_UNKNOWNS * i + 1];
			}
		}
		// Matrix.printMatrix(mainMatrix);
	}

	public double[] getUpdatedBVector(EquationParameters parameters) {
		// TODO implement this
		double[] updatedB = new double[b.length];
		for (Iterator<AbstractElement1D> iterator = elementsList.iterator(); iterator.hasNext();) {
			AbstractElement1D element = iterator.next();
			// element.generateAll();
			ElementMatrix data = element.generateUpdatedVector(parameters);
			List<Integer> pos = data.getNodeList();
			for (int i = 0; i < pos.size(); i++) {
				updatedB[pos.get(i) * NUMBER_OF_UNKNOWNS] += data.vector[NUMBER_OF_UNKNOWNS * i];
				updatedB[pos.get(i) * NUMBER_OF_UNKNOWNS + 1] += data.vector[NUMBER_OF_UNKNOWNS * i + 1];
			}
		}
		return updatedB;
	}

	public void setBothDirihlet(double[] left, double[] right) {
		for (int i = 0; i < NUMBER_OF_UNKNOWNS; i++) {
			mainMatrix.getRows().get(i).getValues().clear();
			mainMatrix.getRows().get(i).getIndexes().clear();
			mainMatrix.setValue(i, i, 1.0);
			b[i] = left[i];

			int nnn = mainMatrix.getN();
			mainMatrix.getRows().get(nnn - i - 1).getValues().clear();
			mainMatrix.getRows().get(nnn - i - 1).getIndexes().clear();
			mainMatrix.setValue(nnn - i - 1, nnn - i - 1, 1.0);
			b[nnn - i - 1] = right[1 - i];
		}

	}

	public static void setBothDirihlet(HashMatrixLUStyle matrix, double[] b, double[] left, double[] right) {
		for (int i = 0; i < NUMBER_OF_UNKNOWNS; i++) {
			matrix.getRows().get(i).getValues().clear();
			matrix.getRows().get(i).getIndexes().clear();
			matrix.setValue(i, i, 1.0);
			b[i] = left[i];

			int nnn = matrix.getN();
			matrix.getRows().get(nnn - i - 1).getValues().clear();
			matrix.getRows().get(nnn - i - 1).getIndexes().clear();
			matrix.setValue(nnn - i - 1, nnn - i - 1, 1.0);
			b[nnn - i - 1] = right[1 - i];
		}
	}

	public double[] solve() {
		LUFactorization lu = new LUFactorization(this.mainMatrix);
		result = lu.calculate(b);
		return result;
	}

	public double[] getB() {
		return b;
	}

	public HashMatrixLUStyle getMainMatrix() {
		return mainMatrix;
	}

	public ArrayList<AbstractElement1D> getElementsList() {
		return elementsList;
	}

	/**
	 * @param args
	 */
	public static void main(String[] args) {
		final IFunction1D f1 = new ConstantFunction1D(0.01);
		final IFunction1D f2 = new ConstantFunction1D(0);
		double v = 0.35;
		double e = 1000;
		double h = 0.01;

		final EquationParameters eqP = new EquationParameters(5.0 / 6.0, e, h, v, f1, f2);
		double a = 0.0;
		double b = 1.0;
		int n = 15;
		double hStep = (b - a) / (double) n;
		final ArrayList<AbstractElement1D> elements = new ArrayList<AbstractElement1D>();
		int i = 0;
		for (i = 0; i < n; i++) {
			// elements.add(new FEElementLinear(eqP, i, a + i * hStep, a + (i + 1) * hStep));
			elements.add(new BubbleElement(eqP, i, a + i * hStep, a + (i + 1) * hStep));
		}

		CylindricMembraneSolver s = new CylindricMembraneSolver(elements);
		s.generate();
		s.setBothDirihlet(new double[] { 0.0, 0.0 }, new double[] { 0.0, 0.0 });
		final double[] res = s.solve();
		// for (int j = 0; j < res.length; j++) {
		// System.out.println(res[j]);
		// }

		XYSeries s1 = new XYSeries("w", true, true);
		XYSeries s2 = new XYSeries("y", true, true);

		// for (int j = 0; j < elements.size() + 1; j++) {
		// s1.add(a + j * hStep, res[j * 3]);
		// s2.add(a + j * hStep, res[j * 3 + 1]);
		// s3.add(a + j * hStep, res[j * 3 + 2]);
		// }

		int NNN = 10;
		AbstractElement1D el = null;
		for (int j = 0; j < elements.size(); j++) {
			el = elements.get(j);
			double hhh = hStep / NNN;
			for (int k = 0; k < NNN; k++) {
				double[] tmp = el.calculate(a + hStep * j + k * hhh, res);
				s1.add(a + hStep * j + k * hhh, tmp[0]);
				s2.add(a + hStep * j + k * hhh, tmp[1]);
			}
		}
		double[] tmp = el.calculate(b, res);
		s1.add(b, tmp[0]);
		s2.add(b, tmp[1]);

		XYSeries s1an = new XYSeries("w an", true, true);
		XYSeries s2an = new XYSeries("y an", true, true);

		final IFunction1D analitW = new IFunction1D() {

			@Override
			public double calculateDerivaty(double x) {
				double p = f1.calculate(x);
				double d = eqP.getD();
				double g = eqP.getG();
				return p / 2.0 / g + 2.0 * (p / 24.0 / d - p / 2.0 / g) * x - p / 4.0 / d * x * x + p / 6.0 / d * x * x * x;
			}

			@Override
			public double calculate(double x) {
				double p = f1.calculate(x);
				double d = eqP.getD();
				double g = eqP.getG();
				return p / 2.0 / g * x + (p / 24.0 / d - p / 2.0 / g) * x * x - p / 12.0 / d * x * x * x + p / 24.0 / d * x * x * x * x;
			}
		};

		final IFunction1D analitGama = new IFunction1D() {

			@Override
			public double calculateDerivaty(double x) {
				double p = f1.calculate(x);
				double d = eqP.getD();
				double g = eqP.getG();
				return -p / 12.0 / d + p / 2.0 / d * x - p / 2.0 / d * x * x;
			}

			@Override
			public double calculate(double x) {
				double p = f1.calculate(x);
				double d = eqP.getD();
				double g = eqP.getG();
				return -p / 12.0 / d * x + p / 4.0 / d * x * x - p / 6.0 / d * x * x * x;
			}
		};

		int nAn = 100;
		double hAn = (b - a) / nAn;
		for (int j = 0; j <= nAn; j++) {
			double xAn = a + hAn * j;
			s1an.add(xAn, analitW.calculate(xAn));
			s2an.add(xAn, analitGama.calculate(xAn));
		}

		XYSeriesCollection col = new XYSeriesCollection();
		col.addSeries(s1);
		col.addSeries(s2);
		col.addSeries(s1an);
		col.addSeries(s2an);
		NewOkCancelDialog d = new NewOkCancelDialog(null, false, col);
		d.setLocationRelativeTo(null);
		d.setVisible(true);

		final IFunction1D wResult = new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double result = 0.0;
				for (AbstractElement1D abstractElement1D : elements) {
					result += abstractElement1D.calculate(x, res)[0];
					if (Math.abs(result) > EPS) {

						break;
					}
				}
				return result;
			}

			@Override
			public double calculateDerivaty(double x) {
				double result = 0.0;
				for (AbstractElement1D abstractElement1D : elements) {
					result += abstractElement1D.calculateDerivaty(x, res)[0];
					if (Math.abs(result) > EPS) {
						break;
					}
				}
				return result;
			}
		};

		final IFunction1D gammaOneResult = new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double result = 0.0;
				for (AbstractElement1D abstractElement1D : elements) {
					result += abstractElement1D.calculate(x, res)[1];
					if (Math.abs(result) > EPS) {
						break;
					}
				}
				return result;
			}

			@Override
			public double calculateDerivaty(double x) {
				double result = 0.0;
				for (AbstractElement1D abstractElement1D : elements) {
					result += abstractElement1D.calculateDerivaty(x, res)[1];
					if (Math.abs(result) > EPS) {
						break;
					}
				}
				return result;
			}
		};

		XYSeries s1func = new XYSeries("gama", true, true);
		XYSeries s2func = new XYSeries("w", true, true);
		s1an = new XYSeries("w an", true, true);
		s2an = new XYSeries("y an", true, true);
		for (int j = 0; j <= nAn; j++) {
			double xAn = a + hAn * j;
			s1func.add(xAn, gammaOneResult.calculateDerivaty(xAn));
			s2func.add(xAn, wResult.calculateDerivaty(xAn));
			s1an.add(xAn, analitW.calculateDerivaty(xAn));
			s2an.add(xAn, analitGama.calculateDerivaty(xAn));
		}

		col = new XYSeriesCollection();
		col.addSeries(s1func);
		col.addSeries(s2func);
		col.addSeries(s1an);
		col.addSeries(s2an);
		d = new NewOkCancelDialog(null, false, col);
		d.setLocationRelativeTo(null);
		d.setVisible(true);

		GaussQuadrature gaussQuadrature = new GaussQuadrature(256);

		double value1 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = gammaOneResult.calculateDerivaty(x);
				double gamma = gammaOneResult.calculate(x);
				double w = wResult.calculate(x);
				double wDer = wResult.calculateDerivaty(x);
				return eqP.getG() * (-gammaDer * gamma + gammaDer * wDer);
			}
		});

		double value2 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = gammaOneResult.calculateDerivaty(x);
				double gamma = gammaOneResult.calculate(x);
				double w = wResult.calculate(x);
				double wDer = wResult.calculateDerivaty(x);
				return eqP.getD() * gammaDer * wDer + eqP.getG() * (gamma * w + wDer * w);
			}
		});

		double value1An = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = analitGama.calculateDerivaty(x);
				double gamma = analitGama.calculate(x);
				double w = analitW.calculate(x);
				double wDer = analitW.calculateDerivaty(x);
				return eqP.getG() * (-gammaDer * gamma + gammaDer * wDer);
			}
		});

		double value2An = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = analitGama.calculateDerivaty(x);
				double gamma = analitGama.calculate(x);
				double w = analitW.calculate(x);
				double wDer = analitW.calculateDerivaty(x);
				return eqP.getD() * gammaDer * wDer + eqP.getG() * (gamma * w + wDer * w);
			}
		});

		System.err.println("v1 = " + value1An + "\t" + value1);
		System.err.println("v2 = " + value2An + "\t" + value2);
		System.err.println("norm = " + (value1An + value2An) + "\t" + (value1 + value2));

		double normL2Cal1 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = gammaOneResult.calculateDerivaty(x);
				double gamma = gammaOneResult.calculate(x);
				double w = wResult.calculate(x);
				double wDer = wResult.calculateDerivaty(x);
				return gamma * gamma;
			}
		});

		double normL2Cal2 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = gammaOneResult.calculateDerivaty(x);
				double gamma = gammaOneResult.calculate(x);
				double w = wResult.calculate(x);
				double wDer = wResult.calculateDerivaty(x);
				return w * w;
			}
		});
		
		double normL2An1 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = analitGama.calculateDerivaty(x);
				double gamma = analitGama.calculate(x);
				double w = analitW.calculate(x);
				double wDer = analitW.calculateDerivaty(x);
				return gamma * gamma;
			}
		});

		double normL2An2 = gaussQuadrature.integrate(a, b, new FunctionWithoutDerivatives1D() {

			@Override
			public double calculate(double x) {
				double gammaDer = analitGama.calculateDerivaty(x);
				double gamma = analitGama.calculate(x);
				double w = analitW.calculate(x);
				double wDer = analitW.calculateDerivaty(x);
				return w * w;
			}
		});
		System.err.println("v1 = " + normL2An1 + "\t" + normL2Cal1);
		System.err.println("v2 = " + normL2An2 + "\t" + normL2Cal2);
		System.err.println("normL2 = " + (normL2An1 + normL2An2) + "\t" + (normL2Cal1 + normL2Cal2));
	}
}